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. j—{ ' Summary. The Yakushevich model provides a very simple pictures of DNA 

1-^ . torsion dynamics, yet yields remarkably correct predictions on certain physical 

qH' characteristics of the dynamics. In the standard Yakushevich model, the in- 

teraction between bases of a pair is modelled by a harmonic potential, which 
becomes anharmonic when described in terms of the rotation angles; here we 
^ . substitute to this different types of improved potentials, providing a more physi- 

■rij" ' cal description of the H-bond mediated interactions between the bases. We focus 

^D . in particular on soliton solutions; the Yakushevich model predicts the correct 

^-P ' size of the nonlinear excitations supposed to model the "transcription bubbles" , 

and this is essentially unchanged with the improved potential. Other features 
^^ ' of soliton dynamics, in particular curvature of soliton field configurations and 

f^ . the Peierls-Nabarro barrier, are instead significantly changed. 

"o 

rJO . Introduction 

. . , In the DNA transcription process, a specialized enzyme (RNA-Polymerase) 

K^ ■ binds to a specific site of the DNA double helix and unwinds it in a local 

k> I region of 15-20 bases, thus creating a "transcription bubble"; the RNAP and 

^ • the bubble travel then along the DNA, copying its sequence and producing a 

5t ! RNA-Messenger to be later used to express genes or replicate the local sequence. 

This process requires a very finely tuned coordination of the motion of RNAP 
- and production of the RNA-Messenger - with the dynamics of the DNA double 
chain. In a pioneering paper appeared in 1980, Englander, Kallenbach, Heeger, 
Krumhansl and Litwin [8] suggested that nonlinear excitations in the DNA 
double chain could be instrumental in this process and allow the motion of the 
transcription bubble to occur at near-zero energy cost. 

*Work supported in part by the Italian MIUR under the program COFIN2004, as part of 
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In particular, as the fundamental motion undergone by DNA nucleotides in 
this process is a roto/torsional one, they suggested to model the DNA molecule 
as a double chain of coupled pendulums; the relevant nonlinear excitations would 
then be (topological) solitons pretty much like those well known in the sine- 
Gordon equation. We refer to their work for detail. 

Their suggestion was taken up by a number of author, producing several 
variants of their fundamental model^ (these are reviewed in [25], to which the 
reader is referred for detail and discussion). In particular, a simple model was 
put forward by L.V. Yakushevich [23], see also [24, 25, 26] for further results 
and extensions; this is recalled and discussed in sect.l below. 

The Yakushevich (Y) model has its appeal, and at the same time weakness, 
in its simplicity; it is remarkable that despite such simplicity it predicts correct 
orders of magnitude for a number of physical features, such as timescales of small 
oscillations and - most relevant in the wake of the approach of Englander et al. 
- size of the nonlinear excitations: in fact, this corresponds to the experimental 
size of the transcription bubble. 

At the same time, the Yakushevich model is a first-grade (i.e. one degree 
of freedom per nucleotide [24]), idealized (i.e. DNA is considered as a homo- 
geneous chain [24]), and nearly- linear model: interaction between successive 
nucleotides is via a harmonic potential, and the intrapair interaction (i.e. in- 
teraction between bases in a Watson-Crick pair, and though these between the 
corresponding nucleotides) is also modelled by a potential which is harmonic in 
the distance, and becomes anharmonic when described in terms of the relevant 
degrees of freedom, i.e. rotation angles. 

While the first two features are nearly unavoidable (but see [20]) in models 
to be analyzed in analytical - and not just numerical - terms, near-linearity 
raises more perplexity even from the point of view of a mathematical approach, 
in particular since the motions we are mainly interested in are fully nonlinear. 

In this note we consider a modified Y model, in which the intrapair potential 
is no longer near-linear. Our main motivation in this study was to understand 
how much of the prediction of Y model survives a fully nonlinear modelling 
of the intrapair potential. Note in this respect that the Pcyrard-Bishop (PB) 
model [19] considers a fully nonlinear Morse potential for this interaction. In 
this case we should also consider the directionality of H bonding, which plays 
no role in the "straight" motions considered in the PB model and is not taken 
into account by the Morse-PB potential. 

More radical modifications of the standard Yakushevich model, modifying 
the other features mentioned above, could also be considered; these will be 
studied elsewhere [3]. 

As the interaction between bases in a pair is mediated by H bonds, wc argue 
that a more realistic approximation is that of a potential describing H-bond 
interactions. This could be a Morse potential as in the PB model, better if 

^In a related but diflferent direetion, other authors also studied models for radial stretching 
motions of the DNA molecule, which are related in particular to DNA denaturation. The 
main line of research in this direction followed the formulation of the Peyrard-Bishop model 
[19] and extensions thereof [1, 2, 18]. In this note wc will focus on roto/torsional dynamics. 



modified by taking into account directionality effects; or even - as we will do in 
the first instance - given that H bonds results from essentially dipolar forces, a 
dipole-dipole interaction. This is not entirely correct, since as bases rotate the 
charge distribution is also modified, but takes into account several features of 
H bonding, in particular its directionality, and of the motion geometry, which 
are not considered in the standard Y model. In a way, the dipole-dipole poten- 
tial tests an approximation which is orthogonal to the one represented by the 
Yakushevich potential. 

It should be mentioned that this work continues investigation described in 
[13], where we considered how the nonlinear dynamics of the Y model is af- 
fected by dropping the contact approximation £o = 0; as in the present case, we 
have shown there that albeit the linear dynamics of the Y model is strongly af- 
fected by the approximations this model assumes, the physical features of fully 
nonlinear excitations are remarkably robust i.e. survive when we drop those 
approximations (that oi £o = in [13], that of a harmonic intrapair potential 
here). The same does not hold when one touches other features of the model, 
as will be shown in forthcoming work [3, 20]. 

In the following, after describing the standard Y model, including also "he- 
hcoidal" terms [5, 9] (sect.l), we discuss in some detail the case of general 
Y-like models and in particular the reduced equations describing soliton solu- 
tions (sect. 2). These equations depend on the intrapair potential V modelling 
H-bond interactions; we also focus on the special case of soliton solutions of 
topological indices (1,0) and (0,1), see below for definitions. In sections 3-5 we 
consider different cases for the intrapair potential: the dipole-dipole potential 
in sect. 3, the Morse-PB potential in sect. 4, and a Morse potential modified by 
introducing a directional term in sect. 5. In all cases, we use the detailed discus- 
sion of sect. 2 and tackle directly the effective equations obtained there; these 
provide a description of the special soliton equations which can be compared - 
as we do - with those provided by the standard Y model. All this discussion 
is conducted on continuous models, but the DNA chain is intrinsically discrete; 
in sect. 6 we evaluate numerically a relevant parameter providing some measure 
of the relevance of discreteness effects, i.e. the height of the Peierls-Nabarro 
barrier. The final sect. 7 is devoted to a discussion of the results obtained and 
to drawing some conclusions. 
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Figure 1: Left: A base pair in the Yakushevich model. Nucleotides are repre- 
sented as disks of radius r; here A = r-|-£o/2 is the distance between the centers 
of disks and the axis of the double helix. Right: Bases which are half-pitch of 
the helix away and on opposite chains are actually near enough in space, and 
Bernal-Fowler filaments [6] joining them can form; they mediate the "helicoidal 
interactions" described by (1.3). 

1 Yakushevich model 



The Yakushevich model (Y model) is an "ideal" DNA model [25]. That is, 
the DNA double chain is considered to be of infinite length, and all bases are 
considered as being equal. 

In the Yakushevich model, each nucleotide is considered as a single unit, 
and its state described in terms of an angular variable. That is, each nucleotide 
is seen as a rigid disk of radius r which can rotate around an axis and has a 
moment of inertia / around this axis; the rotation of the nucleotide at site i € Z 
on the chain a (a = 1, 2) is described by an angle i? 



(a) 



to a specific base pair, we write simply i3i, ^2 for iJ^- ,'d\ 



When we are referring 
The axes of rotation 
of the two disks lie a distance A> r away. 

The model is better described by a Lagrangian L = T — U . In this, the 
kinetic energy is obviously given by 



T-HH^^^ 



q(a)\2 



(1.1) 



The potential energy U is the sum of two terms, or three if we consider the 
so-called "helicoidal" version of the model [5, 9, 25]. These are the stacking 
{Us) term, the pairing one (C/p), and maybe the helicoidal one Uh] the latter is 
relevant only in small amplitude dynamics and in the dispersion relations. Thus 



U = Us + Uh + Uj, 



The first two terms correspond to harmonic potentials^ and depend on di- 
mensional coupling constants Kg and Kh- The first, describing backbone torsion 
and stacking interactions between bases at nearby sites on the same chain, is 
given by 

a i 

The helicoidal term describes interactions between bases which are nearby 
in space due to helical geometry, and mediated in real DNA by Bernal-Fowler 
filaments (chains of hydrogen bonded water molecules [6]), which join a base and 
the one on the other chain half-pitch of the helix away, see fig. lb; in B-DNA, 
this means a distance of ft, = 5 sites along the DNA double chain. These are 
very weak interactions, and are also modelled by a harmonic potential: 

U>^ = ^HHi/t^n-^r^^ ■ (1-3) 

a i 

This term has no effect on fully nonlinear dynamics [14], but introduces qual- 
itative changes in small amplitude dynamics and dispersion relations; we will 
thus drop it when considering fully nonlinear excitations. 

As for the Up term, it describes the nonlinear interactions between bases in 
a pair, i.e. 

Up ^Y. ^(^^^^f^) (1-4) 

i 

with V the intrapair potential. 

We will call the models of the class described so far, with general V , Yaku- 
shevich -like models; we denote the original Yakushcvich model (with V given 
by (1.5) below) as the "standard" Y-modcl. 

In the standard Y model it is assumed that 

V = {l/2)Kp{l-lo? , (1.5) 



where £ is the distance between the atoms bridged by the H bonds, and Iq 
is this distance in the equilibrium configuration described by i?^" = 0. Note 
r + £o/2 = A. 

This harmonic potential becomes anharmonic in terms of the angles "ff^ ■ 
Indeed, with simple algebra (and writing, for ease of notation, -da for ^^° ) we 
have^ e = {£l + Mor + 6r^) - 2r(4 + 2r)[cosz?i + cos??2] + ^r^ cos(i?i - ■&2)- 

In the usual treatments of the Y model, one considers the approximation 
^0 = 0, also called "contact approximation" (a discussion of the usual Y model 
beyond the contact approximation is given in [13]). With this, r ~ A and 

(£-eof = £^ = 2r2[3-2(cos(??i) + cos(^2)) + cos(??i-??2)] ; 



^We stress that the choice of harmonic potentials for these term, and in particular for the 
longitudinal interactions, is common to nearly all DNA models [18, 25]. 

^It should be stressed that — as remarked by Gonzalez and Martin-Landrove [15] — if we 
expand this in a power series for iSi , i?2 around i?i = i?2 = 0, and assuming Ig j^ 0, the 
quadratic term vanishes and (£ - £0)^ « [r2/(4£2)][r(i9i + 1^2)^ + ioi'Sf + 1^I)]^- 



expanding around the equilibrium we get £^ « r^(?9i + ??2)^- Thus, 



V = Ky 



cos(i?f^-i?f^)-2 (cosi?f^+cos^f') +3 . (1.6) 



The intrapair potential term Up reads from this and (1.4). This completes the 
description of the Yakushevich model. ^. 

In studying the Y model it is convenient to pass to variables 

^n = {€' + €^) /2 , Xn = {€^ - €^) /2 (1.5) 

(thus dn — '4'n + Xn and "(?„ — tpn ^ Xn) , Siiid introduce 

a = Ks/I , b = Kh/I , c = 2r^Kp/I) . (1.6) 

The Euler-Lagrange equations stemming from the Yakushevich Lagrangian 
are 

^„ = a {xpn+l - 2lpn + 1pn-l) + b{lpri+h - ^i^n + i^n-h) + 

- csintp cosx , /^ r,N 

Xn = a(X«+l -2x„+Xn-l) - b{Xn+h + '2Xn + Xn-h) + 

— c sin X (cos "0 — cos x) ■ 

Linearizing these about the equilibrium ipn ^ Xn = (^ and writing tpn = 
fqu,exp[i{qn + LL)t)], Xn =■ Sqw cxp[z(gn + a;t)], q = kd, yields the dispersion 
relations 

a;2 = 2a(l-cos(fc(5)) + 2b{l - cos{k6h)) + c; 

col = 2a (1 - cos(fc(5)) + 2b{l + cos (kSh)) . ^ '' 

As for the dimensional parameters appearing in the model, we adopt those 
suggested in [12], i.e. 

r = A = 5.5A, (5 = 3.4 A; / = 3 • lO^-'^^cm^ • g ; 
Ks == 0.13eV/rad^ , Kpr"^ = 0.025eV/rad^ , Kh == 0.009eV/rad^ . ^ ' ' 

The choice for Kp corresponds to an energy of O.leV for the breaking of 
the intrapair coupling (this is AE ~ F(7r,7r) — y(0,0)); this corresponds to 2.5 
times the energy needed for breaking a H bond, i.e. is intermediate between the 
A-T and the G-C cases. 

The full equations admit travelling wave solutions; these are better dis- 
cussed in the continuum approximation in which the arrays {di ,'i?^ }, and 
hence {tpi{t)} and {xi{t)}j are promoted to fields ^a{x,t), <f){x,t), £,{x,t). The 
finite energy condition require that each "daix, t) goes to a multiple of 27r (with 
vanishing x derivative) for \x\ -^ oo; hence finite energy travelling wave solutions 

■^In the following we will consider a different expression for V , which we believe is physically 
better justified than (1.5), which will also dispense us with the prescription io = 0. We will 
first describe other features of the standard Yakushevich model 



are indexed by two integer rotation numbers {m,n). It turns out that there is 
a maximum speed vq = ^/aS for these. 

The fundamental soUton solutions can be easily determined analytically. 
They are given, respectively, by 

Lp{z) = 4arctan[e''^] , ^ = 0; 

^(z) ~ 2 arccos (— kz/\/1 + k'^z^'j , (p = . 



Here we have written k := y^c/{aS^ — Iv^): for v = and with the above choices 
for the parameters, this yields 



K = 5-^ ^2Kpr^/K, w 0.625"^ « 0.18 A ^ . (1.11) 

We refer e.g. to [14, 25] for further detail. 

2 General Y-like models 

As mentioned above, a general Y-like model is obtained when we substitute (1.5) 
with a general intrapair potential V^-di, ^2) having a minimum in t?i = t?2 = 
(and only there, up to 27r periodicity), and symmetric under the exchange of 
the two chains, y(i?i,'i?2) = V'(^2,^i). 

The analysis conducted for the standard Yakushevich model can be to a 
certain extent generalized to this wider class of models. From now on, we will 
use the term "Yakushevich model" in the generalized sense. 

2.1 Equations of motion 

The equations of motion - i.e. the Euler-Lagrange equations arising from the 
modified Lagrangian - are 



(2.1) 



i&:^ = Ks {^\:u ^^'n' + €\ + kh <ri, - 2^^^^ + k-h 



We now pass to the ip^x variables, and write Va for V expressed in terms of 
these to avoid any misunderstanding, i.e. 14 (V', x) = ^("0 + X) V" ~ x)- We get 
(note the sign differences in the Kh term and the 1/2 factor in the Va terms) 

lipn = Ks {tpn+l - 2V'« + V'«-l) + Kh {lpn+h - '^Ipn + 1pn-h) + 

-(1/2)(9K(^,X)/9V^) , (2 2) 

IXn = Ks iXn+l-'^Xn + Xn-l) - Kh [Xn+h + '^Xn + Xn-h) + 

-il/2){dVai^,x)/dx) ■ 

2.2 Small oscillations 

In order to analyze small amplitude dynamics around the 'd^' = equilib- 
rium, we consider the quadratic approximation Vq for V . Due to the symme- 
try property required of V, see above, this will be decomposable as a sym- 
metric and an antisymmetric part for i?i ^ ?92; hence we can write Vq = 



c+{'di + 1)2)'^ + c_(?9i — ?92)^ or, in terms of the tp, x variables, 

K,o(V',x) = (l/2)[(8c+)^2 ^ (8c_)x^ (2.3) 

The dispersion relations are obtained setting "tpnit) = fqu:Cxp[i{qn + cut)], 
ipn{t) = Qquj &xp[i{qn + Lut)], and working at the linear level. In this way, (2.2) 
yields 

w| == 2a[l-cos(<7)] +2fe[l-cos((7/i)]+c^ 

col ^ 2a[l-cos(q)] + 26[l + cos(q/i)]+Cx ^ ' ^ 

for the -0 and the x equations respectively; in these we have used a and b as 
defined in (1.6), and introduced the constants 

c^ := 8c+/I , cx := 8c_// . (2.5) 

Note that the tp branch is acoustical (i.e. limg^o^('?) = 0) or optical (i.e. 
\imq^oOj{q) ^ 0) according to the vanishing or otherwise of the constant c^, 
while for the x branch this depends on the constant Ab + c^. 

2.3 Continuum approximation 

We will now consider the full nonlinear dynamics (2.2). In the study of the 
standard Yakushevich model, the "helicoidal" term turns out to be irrelevant in 
the fully nonlinear regime; on the other hand, they make difficult a discussion 
of soliton solutions. We will therefore drop these (i.e. set Kh = 0) from now on. 
Having set Kh — 0, eqs. (2.2) reduce to 

Ii>n = K, (^„+i - 27A„ + 7A„_i) - (1/2) (9K(V',X)/5V') , (0 P,\ 

IXn = Ks (x„+i-2x„ + Xn-i) - (l/2)(aK(V',x)/9x) • ^ ' ' 

In discussing these, it is convenient to promote the arrays {ipnit)} and 
{Xn{t)} {n G Z) to fields ^{x,t) and E{x,t), the correspondence being given 

by 

Mt) ~ $(n(5,i), xn{t) o^ E{n5,t) . (2.7) 

Here 6 represents the spacing between successive sites; in B-DNA, S ~ 3.4 A. 
With this, and V = (1/2)K, (2.6) reads 

/$„(x, t) = Ks mx + S,t)- 2$(x, t) + $(x - S, t)) - (9F($, S)/a$) {x, t) , 
/S„(x, t) = Ks (S(x + (5, t) - 2S(a;, t) + E{x - S, t)) - (ai>($, S)/as) {x, t) . 

If now we assume that $ and .=, vary slowly in space compared with the 
length scale set by lattice spacing, we can write 



$(x ±S,t)^ $(x, t) ± 5$,(x, t) + (5V2)$,,(x, t) , 
E{x ±S,t)^ S(x, t) ± (5S,(x, t) + ((5V2)S,,(x, t) . 



2/o^^ ^. .^ (2-9) 



Inserting these into (2.8), we finally obtain the field equation for the im- 
proved Yakushevich model in the continuum approximation: 



IEu{x,t) = K,S^E,,{x,t) - (dV{<P,E)/dE^{x,t) . 



(2.10) 



We can from now on omit to write the point (x, t) at which functions should be 
evaluated, as all of them refer to the same point. 

The (2.10) being PDEs, they should be supplemented with a side condition 
determining the function space to which the acceptable solutions belong. The 
physically natural condition is that of finite energy; that is, we should require 
that the integral 



i(/($? + s?) + if.(<i>2 + s2)) + i/($,s) 



da; (2.11) 



is finite; note that if this condition is satisfied at t = 0, it will be so for any t. 

It should also be mentioned that considering a continuum version of the 
model leads to disappearance of the Peierls-Nabarro barrier [18] for soliton mo- 
tion. 

2.4 Travelling wave solutions 

Next we focus on travelling wave solution for (2.10). That is, we restrict (2.10) 
to a space of functions 

$(x,t) ^ ip{x-vt), E{x,t) = -qix-vt) (2.12) 

(we should further restrict this in order to take into account the finite energy 
condition). We will also write simply z :— x ~ vt, and introduce the parameter 

fi := Iv^ - KsS'^ . (2.13) 

With this, (2.10) reduce to two ODEs, i.e. 

^" = - {dW{ip,r^)/dip) , 

r;" = - {dW{ip,J])/dTj) , ^ > 

where we have of course defined 

W{^,rj) := ^i-'V{^,v). 

Tracing back the definition of W in terms of the original intrapair potential V, 
we obtain 

Wiv,il) = ^V{ip + v,^-r^) . (2.15) 

That is, travelling wave solutions are described by the motion of a point 
particle of unit mass in the potential W, with z playing the role of time for this 
motion. Note that /x could be negative, which will actually be the relevant case. 
The conservation of energy reads 

Ui^r + iv'f) + W{^,r^) = E. (2.16) 



2.5 Soliton solutions 

Let us come to the finite energy condition (2.11). In terms of our functions 
(p{z), ri{z), we require that 

Hm tp'{z) = lim r]'{z) = ; (2.17') 

z — ^±oo z — *±oo 

moreover, the functions (p and rj themselves should go to a point of minimum 
for the potential V. 

Note that if /x > 0, the minima of V are the same as the minima of W; but if 
/U < 0, then minima of V are the same as the maxima of W. As one cannot have 
nontrivial motions which reach asymptotically in time (that is, in z) a minimum 
of the effective potential, in order to have travelling wave solutions satisfying 
(6.6') and going to minima of V ior z ^ ±oo we need ji < 0; we assume this 
from now on, and write 

^i = -l/r^ . (2.18) 

In turn, /i < implies that there is a maximum speed for travelling waves: 



|i;| < {VKJi)S ■■= vma. . (2.19) 

The minima of V, i.e. the maxima of W , are for t?^ = 27rcri, i.e. ip ~ 2m:, 
rj = 2TO7r with m,n half- integers, n = (si + S2)/2, m ~ (si — S2)/2; thus the 
finite energy condition requires 



limz^ioo fiz) = 2n±n , limz_»±oo rji^) = 2m±7r 
n± = |(sf + sf) , m± = i(sj^ - sf) ■ 



(2.17") 



The integers <Ti := sj' — s^ are rotation numbers, counting the number of 
complete turns made by the angles i?i in between z = — oo and z = +oo; 
they thus represent a topological index and separate functions satisfying the 
finite energy condition into distinct topological sectors. The same applies to 

n := n+ — n_ and m := m^ — 77i_ 

(we can always take n_ = m^ = by suitably choosing the origin for the angles 
if and ry). The solutions with nonzero (rn,n) will correspond to topological 
solitons. As the problem admits a variational formulation, we are sure there are 
solutions in each topological sector [7]. 

In terms of the dynamical system describing the evolution in "time" Z in the 
potential W, they represent heterochnic solutions connecting the point (0, 0) at 
z = — oo with the point (27rn, 2Tnn) at z ~ +oo. It is thus no surprise that the 
analytic determination of such solutions is in general impossible. 

2.6 Special soliton solutions 

The soliton solutions with indices (1,0) (i.e. cti = (T2 = 1) and (0,1) (i.e. 
(Ti = 1 = — (T2) are special and can be determined explicitly. 
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These are special in that they make vary only one of the two fields. That is, 
the (1,0) solution will have r]{z) = 0, and the (0, 1) solution will have (p{z) = 0. 
Thus, they correspond to one-dimensional motions in the effective potential 
W{(f,i]), and as such they can be exactly integrated, as we now discuss. 

By adding a constant we can always set V{0, 0) = and hence W{0, 0) — 0, 
see (2.15), so we assume this to be the case. 

The (1,0) soliton corresponds to 77 = 0, so we define 

P{ip) :- W{^,0) , (2.20) 

and rewrite the conservation of energy (2.16) as 

We are interested in the motion with E = Eq = W{Q, 0); we assumed this to be 

£^0 = 0. 

Eq.(2.21) yields then the separable equation 



^ ^ V-2P(^), (2.22') 

or equivalently 

-^^L= = dz . (2.22") 

Introducing the reduced effective potential (the second equality follows from 
(2.15)) 

P{ip) :- 2^iP{^) - V{^,^) (2.23) 

and recalling (2.18), we get immediately the equation for the (1,0) soliton: 

dtp 



VpW) 

Denoting by / the integral of the l.h.s. 



Tdz. (2.24) 



/(V»:^/^, (2.25) 

we get hence 

/(^) = r(z - zo) , V = rM^^ - ^o)] • (2.26) 

The same approach allows to obtain the (0,1) soliton, for which (p = Q. In 
this case we define 

Q{r^) ~ VK(0,r7), (2.27) 

and rewrite the conservation of energy (2.16) as 

i(g)\0(,) = i.. (2.28) 
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We again consider E = Eq = 0: eq.(2.28) yields then 

"^2= = dz , (2.29) 

and introducing, see again (2.15) for the second equahty, 

Q{r,) := 2^iQ{f1) = V{r,,-r,) , (2.30) 

we get immediately the equation for the (0,1) soliton: 

dip 



Denoting by g the integral of the l.h.s., 



Tdz. (2.31) 



we get hence 

g(?7) = t(z - zo) , 77 = 5"^[t(z-zo)] • (2.33) 

It may be worth recalling, as a final remark in our general discussion, that 
in view of (2.23) and (2.30), the equations (2.22) and (2.31) obeyed by the (1,0) 
and the (0,1) solitons are also written in terms of the original intrapair potential 
V(t?i,x?2)as 

d^/dC - Vn^ (2.34) 

for the (1,0) soliton, and 



dr;/dC = VV{r],-v) (2.35) 

for the (0,1) soliton. Here we introduced the rescaled variable 



C := T z = z/^\ = z/(S^Ks - I{vl5f) ; (2.36) 

the latter equality follows from (2.13). 

3 Improving the intrapair potential I. 
Dipole-dipole interaction 

The intrapair potential V is due to the H bonds between complementary bases. 
As H bonding is ultimately a dipolar interaction, on a physical basis it ap- 
pears preferable to consider V as, & dipole-dipole potential rather than simply 
a harmonic potential depending on variations of the euclidean distance of the 
involved atoms. 

Our first improved model will therefore be described as follows: on the edge 
of the disks representing nucleotides in the Yakushevich model sits dipole of fixed 
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strength and whose orientation follows the bases orientation. For 191=192= 0, 
the dipole have the same orientation in space; as the angles change, they at the 
same time acquire discordant orientation and vary their mutual distance. The 
potential corresponding to given angles can be easily computed by elementary 
trigonometry. 

3.1 Intrapair potential 

We consider an intrapair potential due to a dipole-dipolc interaction; this is 
obtained by expanding the full dipole-dipolc interaction potential - obtained 
proceeding as described above ~ in terms of the dipole momentum and keeping 
only second order terms. Proceeding in this way we obtain, with 7 a dimensional 
parameter and Cq an arbitrary constant 

V = Co + 7 [14Ar(cosi9i +COS192) + 2Ar[cos(i9i -2i92) +cos(i92 -2i9i)] + 
-4(^2 + 2r^) cos(i9i - i92) - 12^^ cos(i9i + i92) + 
-r2 cos(2(i9i - i92)) - 7r^] x 

X 2V2(2A2^r2(l + cos(i9i -i92)) -2Ar(cosi9i +cosi92))''^''^j 

(3.1) 
We are of course assuming that A ^ r, i.e. we are not adopting the Yakushevich 
approximation £0 = 0. 

The above will be our choice for the intrapair potential, hence we will have 

Up = 5]n^r\i9fV (3.2) 

i 

As for the other potential terms, we retain Yakushcvich's expressions. 
In terms of the ip and x variables defined above, V reads 

Va = Co - 7 [7r2 + 12^42 cos(2i/;) + 4(^2 ^2r2)cos(2x) + r2cos(4x)+ 
— IGArcosi/icosx (1 + cos^ x)] x 

p yo"| — 1 

X 16 (A^ + r^ cos^ X ^ 2j4r cos tp cos x) 

(3.3) 

The potential is plotted in fig. 2. Note that it has the properties we should 

expect from a potential representing H bonding: it is strongly directional, and 

the atoms are essentially free when their position is not near enough to the 

equilibrium position. 

The potential V (or Va) depends on three dimensional parameters, i.e. the 
coupling constant 7 and the two distances A and r. It may be convenient - 
in particular, in order to compare the results for this model with those for the 
standard Yakushevich model - to write it in terms of 7, A and the adimensional 
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Figure 2: Left: The potential V{'di,-d2) as given by (3.1). Here we set 7 = 1 and 
use the geometrical parameters A = 5.5 A, r = 4.0 A. Darker areas correspond 
to lower potential energy. Right: The intrapair potential Vaiip, x) in terms of ip 
and X (parameters as above). Note that beside the 27r periodicity in ip and in x, 
expressed by Va{ip, x) — ^a(V' + 27r, x) = Va{ip, x + 27r), there is also a periodicity 
for simultaneous tt translations in ip and Xj i-e. V^('0iX) — ^a(V' + ''''i X + ""); 
this follows immediately from (1.5). 



constant A = r/A < 1. In this way, (3.1) reads 

V ^ Co + 7 [2A[7(cosi9i +COS192) +cos(i9i -2i92) + cos(i92 -2i?i)] + 
-4(1 + 2A2) cos(i9i - ^2) - 12 cos(i9i + -32) + 
-A2[7 + cos(2(i?i-i?2))]] X 

X 2V2A3(2 + A2(l + cos(?9i-'i92)) -2A(cos??i+cos^2))^^^j 

(3.4) 
We also note that in order to have 1^(0, 0) = 0, the additive constant cq 
should be chosen as 

CO = -j/[A{l-X)f . (3.5) 



3.2 Model parameters 

Our model depends on several parameters. Some of these ~ i.e. A, Kg and Kh 
are common to the standard Yakushcvich model, and we will retain the values 
given in sect.l for these. 

Note that now r is no longer coinciding with A] as the length of H bonds is 
€0 — 3A and 2r — 2A — io, we take r = 4.0A; that is, we choose 



= 4.0A , £0 = 3.OA , A = 5.5A ; A = {r/A) ~ 0.7272. 



(3.6) 



In order to select a value for 7, we note that the equilibrium corresponds to 
di = = i?2, while the "maximally open" state of the base pair is obtained for 
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^1 = Ti = ^2- Thus y(7r, tt) — F(0, 0) should be of the order of the energy A£' 
needed to break the base pair; with the choice aheady accepted for the standard 
Yakushevich model, we should thus have l^(7r,7r) — F(0,0) = AE ~ O.leV. 
On the other hand, (3.1) yields 



F(0,0) = Co - l/{A-rf , y(^,7r) = cq - l/{A + rf ; 



(3.7) 



hence 

F(7r,7r)-F(0,0) = 7 



1 



1 



7 



1 



1 



{A-rY {A + rfJ A"^ \{l - Xf (1 + A)3 

thus 7 should be chosen according to 

{A^ - r2) 



7 = (Ai?) 



.2^.2, = i^E)A 



3 (l-A^) 



2\3 



2r(3yl2+r2) '~^ ' " 2X{3 + X)' 

The relevant parameter for V is 7/4^'^, and by (3.9) we have 

(1-A2)3 



(3.8) 
(3.9) 



7 
4^3 



(Ai^) 



8A(3 + A) 
With the values given above for AE, A and r, this yields 

7w0.34A3eV; ^ « 5.1 • 10"'' eV . 



(3.10) 



(3.11) 



3.3 Special soliton solutions 

As discussed in sect. 2, soliton solutions are obtained as motions in the effective 
potential W{(f,ri), see (2.14); in the present case W is given - as follows from 
(2.15) - by 

W{^,rj) = -co/(2|m|) + 

+ [7/(2|/z|v43)] [7A2 + 12 cos(2i/;) + 4(1 + 2A2) cos(2x)+ 
+A2 cos(4x) -- 16A cos tp cos x (l + cos^ x)] x 

X 16 (1 + A^ cos^ X — 2Acos?/jcosx) 

(3.12) 
It is convenient, for the computations to follow, to set cq as given in (3.5). 

3.3.1 The (1,0) soliton 

The (1,0) soliton is obtained, as discussed in sect. 2, by considering P{(f) := 
V{(f, tf) and solving (2.24). In the present case it results 



^<'" - Ta' 



1 + 3cos(2(p) - 8Acos(v3) + 4A2 



(1-A)3 (l-2Acos^ + A2)5/2 

where we have also used (3.5); we have then to solve 



(3.13) 
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Figure 3: On the left: The potential Va{ip,x) as a function of ip for x — 
(solid line) and as a function of x for ijj = (dashed line); we have used the 
values A = 5.5 A, r = 4 Aand set /i = —1. On the right: The effective 
potentials P(ip) = W(lp,0) (solid) and Q(rj) = T4^(0,?7) (dashed) entering in the 
determination of the (1,0) and the (0,1) solitons. The additive constant cq has 
been chosen here so to have V{0, 0) = 0, and energy is measured in units such 
that y(7r, 7r) - ^(0, 0) = A^ = 1. 



= T^Pi^p), 



(3.14) 



We were not able to perform the integration (2.24) analytically. On the 
other hand the equivalent equation 

dip 

dz 

can be integrated numerically. The result of this integration is shown in fig. 4a, 
where it is also compared with the standard Yakushevich (1,0) soliton. 

3.3.2 The (0,1) soliton 

The (0,1) soliton is obtained, as discussed in sect. 2, by considering Q{ri) := 
V{rj, —rj) and solving (2.31). In the present case it results 



Qiv) = 



7 
4^3 



2(1 



S^l]) 



(3.15) 



(1-A)3 (l-2Acosr/ + A2)5/2_ 

where we have also used (3.5). Again we were not able to perform the integration 
(2.31) analytically, but the equivalent equation 

d?7 

dl 

can be integrated numerically. The result of this integration is shown in fig. 4b, 
where it is also compared with the standard Yakushevich (0,1) soliton. 



Voiv) 



(3.16) 



4 Improving the intrapair potential II. 
Morse potential 

The problem of modelling intrapair interactions mediated by H bonds has been 
tackled in all attempts to provide dynamical models of DNA; among these, the 



16 









6 


_/-'7 












5 


/ / 












4 


^ 












/f 












/ 


,.1 








-30 


-20 


-10 




10 


20 


30 









6 











■■''/^^' 












5 














4 














A' 












J- 


/''i 








-30 


-20 


-10 




10 


20 


30 



Figure 4: Comparison of the (1,0) (left) and the (0,1) (right) soUtons for the 
standard Yakushevich model (dotted) and our modified model (solid lines). De- 
spite a certain difference in shape, the size of the (1,0) soliton - corresponding to 
symmetric motions in the two helices - is remarkably similar in the two models. 





Figure 5: Left: The Morse-PB potential (4.1); here Z) = 1, /c = 1 xo = 1, and 
we add a constant so that F(xo) — 0. Right: the resulting potential in terms of 
the rotation angle -d, for dx — -62 (solid) and for 'd\ = —192 (dashed). 

Peyrard-Bishop model (and later improvements by Dauxois [5] and by Barbi, 
Cocco and Peyrard [1, 2]) has a special relevance in view of its success in de- 
scribing DNA denaturation [18]. In this, intrapair interactions are described by 
means of a Morse potential 



V{r) 



-k{-i — ro) 



(4.1) 



where r is the distance between (reference atoms belonging to) bases in the 
given pair. 

Note that V has a minimum at r = tq with V{ro) — —a and limr^oo ^(^) = 
0, so that a represents the binding energy. We also have V"{rQ) = 2ak^. 

Thus, albeit the Peyrard-Bishop (PB) model aims at describing a different 
process undergone by DNA, it is quite natural to consider their modelization of 
the intrapair potential. 

It should be noted that - as pointed out in (4.1) - the Morse-PB potential 
depends only on the distance between H-bonded atoms, and not on their aligne- 
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ment. This is no surprise, since in the PB model and later improvements the 
bases lie on a straight line passing through the axis of the double helix. ^ 

In our case we have to consider rotational motions, and thus the latter 
conditions is not holding any more. This raises the question of modelling angular 
dependence of H bonds in this context; in this section we will consider the 
isotropic Morse-PB potential (4.1), while in the next section we introduce an 
angular dependence. 

Also, for our purposes it is more convenient to set additive constants in V so 
that V^(0, 0) = 0. Doing this, and denoting as in previous section the distance 
between (relevant atoms in) the two bases of the pair as i, and the distance 
equilibrium as Iq, the Morse-PB potential reads 

F(t?i,7?2) = a (e-fcW''i.''^)-^o)-l)' . (4.2) 

Here we have pointed out explicitly the dependence on the angles 1?^. We also 
write, for ease of notation, 

p{^i,^2) :- 1(^1,^2) - io , (4.3) 

so that (4.2) reads 

T/(^i,^2) = a (e-'^pi^^^^^) -lY . (4.4) 

4.1 Soliton solutions 

With our choice of variables, the distance £ is written as 

£ = 2A[l- A(cos -i^i + cos ^2) + (A^S) (1 + cos(^i - ^2))] ^^^ (4.5) 

and of course £0 ~ 2A{1 — A). 
By standard computation. 



p{ip, (p) = 2A (^/l-2\cosifi + \'^ - (1 - A)] , 
Pifli-V) = '2A {yj\ - 2A cos r/ + A2 cos2 77 - (1 - A)] 



(4.6) 



This and (4.4) immediately provide the symmetric and antisymmetric re- 
ductions, P((/?) and Q(r]), of the intrapair potential which determine the (1,0) 
and (0,1) solitons via (2.34) and (2.35); see (2.10) below. 

^Moreover, they are always assumed to move symmetrically. Asymmetrical motions are 
also possible dynamically, but less favorable energetically; moreover the dynamics is naturally 
split into symmetrical and antisymmetrical motion, due to an exchange symmetry. Thus one 
usually focuses on the more relevant sector of symmetrical motions [14]. 
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4.2 Model parameters 

The intrapair Morse-PB potential Vp and hence the effective potential W depend 
on the two parameters a and fc; this is at variance with the Yakushevich potential 
and the dipole-dipolc potential considered in sect. 3: they both depend on a 
single parameter. 

We can therefore impose two physical conditions in order to fix these pa- 
rameters; we will retain the values adopted above for the other parameters. 

As pointed out above, a represents the dissociation energy for the Morse 
interaction, supposed to model the H-bond; hence we should have 

a ^ AE ~ 0.1 eV . (4.7) 

We also pointed out that for £ ~ £o, V ~ ak'^{£ — ^o)^; comparing this with 
(1.5), we require that 

2ak^ = Kp . (4.8) 

Here Kp should be chosen as in (1.10), taking into account our choice for r. 
This yields 



k = ^Kp/{2a) = r-^ ^Kpr^/i2AE) ~ 0.09A ^ . (4.9) 



4.3 Special soliton solutions 

The relevant reduced potentials are now 



\2 



P(ip) = V{ip,ip) = a (cxp [-2AfcA(Vl - 2A cos ip - {1 - X))] - l) 

Qiv) = V{r],-r]) = a (cxp[-2AfcA(l-cos77)])^ . 

(4.10) 

The equations (2.34) and (2.35) governing the special (1,0) and (0,1) soli- 
ton solutions are again impossible to solve analytically, but can be integrated 
numerically. The results of these integration are shown in fig. 6, where we also 
plot the standard Yakushevich solitons of the same topological number for com- 
parison. Quite remarkably, again the overall size is about the same as with the 
standard Yakushevich model. 

5 Improving the intrapair potential III. 
Morse potential with a directional term 

As mentioned in the beginning of the previous section, the Morse-PB potential 
provides an "isotropic" description of the H-bond mediated intrapair interaction. 
On the other hand, as discussed in the Introduction, once we allow rotational 
motions, hence the possibility that the alignment between the atoms intervening 
in the H-bond is broken, we should consider the directional nature of the H- 
bonds. This makes that the bases are essentially free (as far as this interaction 
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Figure 6: Special soliton solutions for Morse intrapair potential (solid lines) 
compared with standard Yakushevich solitons (dotted). Left: the (1,0) soliton; 
right: the (0,1) soliton. See text for parameters values; we plot solutions for 
w = 0. 



is concerned) once they are more than a few degrees away from the equilibrium 
position. 

It should be stressed that in our model the bases cannot come any nearer 
than their equilibrium distance; that is, we do not have to worry about the 
region on the left of the minimum of Vp{r). 

A simple way to modify the Morse-PB potential in order to take into account 
the angular dependence of H-bonds is to introduce a prefactor ^(^i, ^2) which 
is unity at equilibrium and goes rather quickly to zero when at least one of the 
bases rotate away from this. It is convenient to choose a smooth function for 
this A; it should moreover be 27r-periodic in both its argument. 

The natural candidate for a prefactor with these features is a gaussian-like 
function: we choose (no confusion should be possible with the 7 used in sect. 3) 



A(^i,^2) := exp[-7((l-cos^i) + (l-cosi?2))] 



We thus have a modified Morse-PB potential 



V = A{^i,^2) [a((e 



-fcp(l5i,)?2) 



1) 



1 



(5.1) 



(5.2) 



where we have used notation introduced in sect. 4 above. Adding a constant a 
so that V(0, 0) = 0, the explicit expression for this is 

V = a exp[-7(l-cos^i)(l-cos^2)] \{1 ~ exp[-kp{-di,d2)] - if - 1 

(5.3) 

We will take the same values as in the previous sections for the parameters 
of the Morse potential. As for the newly introduced parameter 7, this should be 
chosen so that the H-bond interaction is conveniently reduced when the relevant 
atoms are not properly aligned. A natural criterion is to require A{9o, 0q) = ag 
with somewhat arbitrary values for ^o ^^nd oq. If we choose 9o ~ 7r/4 and 
oo = 1/10, we get 

7 ~ 4 ; (5.4) 
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Figure 7: The Morse directional potential V{di, 1^2) as given by (5.3) with j ~ 4 

(left); and the same in adapted coordinates V',X: i-^- K(V':X) = ^(V' + Xi^^x) 
(right). These can be compared with fig. 2. 

this will be our choice. The resulting potential is plotted in fig. 7. 

Special soliton solutions are obtained as usual by numerically solving equa- 
tions (2.34) and (2.35). The relevant reduced effective potentials are in this 
case 



Qiv) 



1 + g-27(i-cos¥^) ((exp[-2fcAF] - f)2 - 1)] 
f + e-27(i-cosr,) ((exp[-2fcAG] - If - 1)] 



where we have defined, for ease of writing, 



F 



v/l-2Acos^ + A2 - (1 - A) , G = A (1 - cos?7) 



(5.5') 



(5.5" 



Needless to say, our choice for the parameter 7 was to a large extent arbi- 
trary; we should then check that the resulting soliton solutions are not strongly 
dependent on such arbitrary choice. To this purpose, we have integrated the 
equations also using different values for 7, and it resulted that the soliton so- 
lutions arc very little dependent on this. The results of integration with 7 = 1 
and for 7 = 10 are shown, together with the solution obtained for 7 = 4, in 
fig.9. 

6 Discreteness effects 

We have described solitons in the continuum approximation; our original model 
was however a discrete one and as well known, going back to a discretized 
description the soliton propagation results to be hindered by the Peierls-Nabarro 
barrier [18]. 
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Figure 8: The (1,0) soliton (left) and the (0,1) soliton (right) solutions for 
the Morse intrapair potential with directional term (5.3) (solid lines) and, for 
comparison purposes, for the standard Yakushevich model. Here parameters of 
the Morse potential and geometrical parameters are as in section 4, i.e. A = 5.5 
A, r == 4.0 A, A = 0.73, a = 0.1 eV, k = 0.088 A"\ and 7 = 4. 





Figure 9: Dependence of the (1,0) (left) and (0,1) (right) soliton solutions for the 
Morse intrapair potential with directional term on the value of the parameter 7. 
Keeping all other parameters equal (values as in fig. 8 above), we plot solutions 
for 7 = 4 (solid line) together with those for 7=1 (dashed) and for 7 = 10 
(dotted). Solutions are extremely similar. 
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The intrapair potentials considered in this note lead to effective potentials 
which are steeper than the one corresponding to the original Yakushevich po- 
tential near the equilibrium position, but have a nearly flat region at larger 
amplitude of the 0i fields; moreover, they lead to soliton solutions which have a 
much larger curvature at the interface between the nearly quiescent region and 
the "transition" region (in which the field passes from to 27r), and a nearly 
zero curvature in the transition region. These effects lead to opposite variation, 
and so the net effect should be evaluated quantitatively in order to understand 
which contribution is larger. 

We thus want to evaluate the Peierls barrier for the different potentials we 
have considered; we will focus on the (1,0) special soliton solutions. The Peierls 
barrier for a given travelling wave can be evaluated as follows. 

First of all, we can just consider static solutions (/'(z) = f{x), i.e. set v = 0; 
for these the energy of the configuration will not include the kinetic term, and 
reduce (see sect.l) to 

H ^ Us +Up . (6.1) 

Note that the balance between these two contribution is only marginally de- 
pendent on the model we consider: indeed parameters were determined by the 
strength of the stacking and the pairing interactions, as they result from exper- 
iments. 

The field ip{x) should now be replaced by the array of values 

^n ■■= <p{nS) . (6.2) 

The stacking energy U^ is given by (1.2). As in our case ^i{x) = ^2{x) = 
ip{x), the latter reads simply 

Us ^ KsJ2 ('^"+1 - '^»)' ; (6-3) 

n 

note this holds for any Y-like model, and in particular for all the models we 
have considered. 

The pairing energy Up depends on the intrapair potential V , and is given, 
see again sect.l and in particular (1.4), by the sum of the pairing energy of all 
base pairs. In our case. 

Up ^ Y. Vi^n,^n) - Y. ^('^") ■ (6.4) 

n n 

Using (6.2), we can express the energy for the chain configuration in terms 
of the soliton field configuration; indeed (6.3) and (6.4) are also written as 

Us = KsJ2 ["^((^ + 1)'^) - "^M)]" ^ Up = J2 ^bM)] • (6-5) 

n n 

We introduce now a shift parameter s G [0, 1], and compute the energy when 
the configuration is shifted by sS: 

H{s) := Ks ^[(/^((n+l + s)(5)-(^((n + s)(5)]2 + ^ P[</j((n + s)5)] . (6.6) 
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The Peierls barrier is then given by the different between the maximum and the 
minimum over s G [0, 1] of the function 



B{s) := H{s) - H{0) 



(6.7) 



We have computed this for the different models considered here; the results 
are summarized in the table below: 



Intrapair potential 


Peierls barrier 


Standard Y 




4.3 • 10-=^ eV 


Dipole 




2.3 • 10-2 cV 


Morse-PB 




3.2 • 10"^ eV 


Directional Morse (7 = 


= 4) 


1.1 • 10-2 eV 


Directional Morse (7 = 


= 1) 


3.5 • 10--* eV 



It results that considering directional rather than isotropic potentials (that 
is, dipole-dipole rather than the standard Y potential; or directional Morse 
rather than isotropic Morse-PB) leads to a higher Peierls barrier, i.e. isotropic 
potentials make us underestimate its height. 

On the other hand, a Morse-type potential (Morse-PB rather than standard 
Y; or directional Morse rather than dipole-dipole) lead to a lower Peierls barrier. 

7 Discussion and conclusions 

We have considered the Yakushcvich model for DNA torsion dynamics, in a 
version modified by introducing different expressions for the intrapair potential 
V{'&i,'&2); these were: 

(a) a dipole-dipole interaction, which takes into account the essentially dipolar 
nature of the H bonds mediating the intrapair interaction between bases in a 
pair; 

(b) a Morse-type potential, as often used to model H-bonds; 

(c) a Morse-type potential modified by introducing a prefactor which takes into 
account the directional nature of the essentially dipolar interaction at the origin 
of H-bonds. 

We have focused on travelling soliton solutions, in particular the special ones 

- with topological indices (1,0) and (0,1) respectively - in which only one of the 
fields ip and x varies from the equilibrium configuration. 

We found that these solutions vary quite little when we consider these dif- 
ferent expressions for the potential (on the other hand, the Peierls barrier is 
quite different in the different cases considered). In particular, the width of the 
soliton is, provided we choose parameters appearing in the model following the 
same physical criteria, of the same order of magnitude. 

This should not be too surprising: after all, once we fix boundary conditions 

- i.e. the topological indices - the solution results from a balance between 
the stacking interaction (which would favor very smooth transitions from one 
vacuum to the other) and the on-site intrapair potential (which would favor as 
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little bases as possible being away from minima of this potential). Thus, the 
width of the soliton depends essentially on the ration between strength of these 
two interactions. 

The Yakushevich model is concerned with DNA torsion dynamics; according 
to the pioneering proposal of Englander et al. [8] , topological (kink) solitons in 
this dynamics should somehow correspond to the "transcription bubbles" which 
are formed when RNA-Polymerase binds to the DNA double chain; from this 
point of view, the main success of the Yakushevich model was to provide the 
correct size for these nonlinear excitations. 

We found that the soliton solutions have a different shape in the standard 
and in the modified Yakushevich models; but their size is not very different in 
the different models. We have also computed the Peierls barrier for the different 
models, finding a quite large range of variability. 

In conclusion, our study suggests that: 

• (A) The standard Yakushevich model captures to a large extent the essen- 
tial part of DNA dynamics as it can be described by simple models with 
a single degree of freedom per nucleotide. 

• (B) On the other hand, study of the standard Yakushevich model can have 
led to an overestimation of the Peierls-Nabarro barrier to be overcome by 
the topological solitons to move along the discrete DNA double chain. 

Point (B) suggests that solitons with too low speed could actually be unable 
to propagate along the chain, and thus set also a minimum speed - and not just 
a maximum one, see sect. 2 - for such solutions, at variance with the analysis of 
purely continuous solitonic excitations in the Yakushevich model [11]. 

Point (A) shows that it is pointless to attempt an improvement of the stan- 
dard Yakushevich model - e.g. to reconcile the predictions of the model with 
different physical requirements, which leads to contrasting tuning of the param- 
eter values [26] - simply by improving the expressions of the intrapair (pairing) 
potential. 

Albeit we have not considered this point here, the results of this study and 
the qualitative discussion above seem to suggest that the same applies to some 
extent for improved expressions of the stacking potential. In this respect we 
recall, however, that recent work by Saccomandi and Sgura [20] points out new 
phenomena appearing as a nonlinear stacking potential is considered. 

Thus the main result of our investigation is that in order to go over the 
limitation of the standard Yakushevich model we should not try to investigate 
more and more detailed description of the interactions between different parts of 
the DNA molecule within the description provided by the standard Yakushevich 
model, but rather consider a more detailed model, maybe with equally rough 
approximation of the interactions. In the language of [24], we should go to a 
higher level in the hierarchy of DNA models. 

This will be done in a different note [3], where a Y-typc model with more 
degrees of freedom will be considered. On the other hand, the investigation 
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reported in the present paper suggests that the new model can be tackled, 
in the first approximation, considering very simple expressions for interaction 
potentials. 

References 

[1] M. Barbi, S. Cocco and M. Peyrard, "Helicoidal model for DNA opening", 
Phys. Lett. A 253 (1999), 358-369 

[2] M. Barbi, S. Cocco, M. Peyrard and S. Ruffo, "A twist-opening model of 
DNA", J. Biol. Phys. 24 (1999), 97-114 

[3] M. Cadoni, R. De Leo and G. Gaeta, "A composite Y model of DNA torsion 
dynamics" ; forthcoming paper 

[4] C. Calladinc and H. Drew, Understanding DNA, Academic Press (London) 
1992; C. Calladine, H. Drew, B. Luisi and A. Travers, Understanding DNA 
{3^"^ edition). Academic Press (London) 2004 

[5] Th. Dauxois, "Dynamics of breather modes in a nonlinear helicoidal model 
of DNA", Phys. Lett. A 159 (1991), 390-395 

[6] A.S. Davydov, Solitons in Molecular Systems, Kluwer (Dordrecht) 1981 

[7] A. Dubrovin, S.P. Novikov and A. Fomenko, Modern geometry. Springer 
(Berlin) 1984 

[8] S.W. Englander, N.R. Kallenbach, A.J. Heeger, J. A. Krumhansl and A. 
Litwin, "Nature of the open state in long polynucleotide double helices: 
possibility of soliton excitations", PNAS USA 77 (1980), 7222-7226 

[9] G. Gaeta, "On a model of DNA torsion dynamics", Phys. Lett. A 143 
(1990), 227-232 

[10] G. Gaeta, "Solitons in planar and helicoidal Yakushevich model of DNA 
dynamics", Phys. Lett. A 168 (1992), 383-389 

[11] G. Gaeta: "A realistic version of the Y model for DNA dynamics; and 
selection of soliton speed"; Phys. Lett. A 190 (1994), 301-308 

[12] G. Gaeta, "Results and limits of the soliton theory of DNA transcription", 
J. Biol. Phys. 24 (1999), 81-96 

[13] G. Gaeta, "Solitons in the Yakushevich model of DNA beyond the contact 
approximation" , preprint 2006 

[14] G. Gaeta, C. Reiss, M. Peyrard and Th. Dauxois, "Simple models of non- 
linear DNA dynamics", Rivista del Nuovo Cimento 17 (1994) n.4, 1-48 



26 



[15] J. A. Gonzalez and M. Martin-Landrove, "Solitons in a nonlinear DNA 
model", Phys. Lett. A 191 (1994), 409-415 

[16] R. Lavery, A. Lebrun, J.F. AUemand, D. Bensimon and V. Croquette, 
"Structure and mechanics of single biomolecules: experiments and simula- 
tion", J. Phys.: Condens. Matter lA (2002), R383-R414 

[17] M. Peyrard (editor). Nonlinear excitations in biomolecules, (Proceedings of 
a workshop held in Les Houches, 1994), Springer (Berlin) and Les Editions 
de Physique (Paris) 1995 

[18] M. Peyrard, "Nonlinear dynamics and statistical physics of DNA" , Non- 
linearity 17 (2004) R1-R40 

[19] M. Peyrard and A.R. Bishop, "Statistical mechanics of a nonlinear model 
for DNA denaturation" , Phys. Rev. Lett. 62 (1989), 2755-2758 

[20] G. Saccomandi and I. Sgura, "The relevance of nonlinear stacking interac- 
tions in simple models of double-stranded DNA" , preprint 2006 

[21] W. Saenger, Principles of nucleic acid structure, Springer (Berlin) 1984 

[22] S.B. Smith, L. Finzi and C. Bustamante, "Direct mechanical measurements 
of the elasticity of single DNA molecules by using magnetic beads" , Science 
258 (1992), 1122-1126 

[23] L. V. Yakushevich, "Nonlinear DNA dynamics: a new model" , Phys. Lett. 
A 136 (1989), 413-417 

[24] L.V. Yakushevich, "Nonlinear DNA dynamics: hyerarchy of the models", 
Physica D 79 (1994), 77-86 

[25] L.V. Yakushevich, Nonlinear Physics of DNA, Wiley (Chichester) 1998; 
second edition 2004 

[26] L.V. Yakushevich, A.V. Savin and L.I. Manevitch, "Nonhnear dynamics of 
topological solitons in DNA", Phys. Rev. E QQ (2002), 016614 



27 



